function calcMovieKernel(R)

% R: N by M rating matrix, N = #users, M = #movies

N = size(R, 1);     % #users
M = size(R, 2);     % #movies

K_v = zeros(M,M);

% The generator matrix
W = zeros(M,M);

% Compute the generator matrix
for i = 1 : M
    for j = i : M 
        
        % Cosine similarity
        W(i,j) = R(:,i)' * R(:,j) / (norm(R(:,i)) * norm(R(:,j)));
    end
    W(:,i) = W(i,:);
end

%*******************************************************************
% Construct the epsilon-neighborhood graph
thres = 0.3;
for i = 1 : M
    ind = W(i,:) >= thres;
    W(i,:) = 0;
    W(i,ind) = 1;
    W(i,i) = 0; % FIXME: necessary?
    W(:,i) = W(i,:);
end

%*******************************************************************
% Exponential kernel
K_v = expm(0.2 * W);
K_v_inv = inv(K_v);

save('u1_movie_exp_kernel.mat', 'K_v');
save('u1_movie_exp_kernel_inv.mat', 'K_v_inv');

%*******************************************************************
% Diffusion kernel

% beta = 0.005

% % The degree matrix
% Deg = diag(sum(W));
% 
% % Negation of graph laplacian
% L_v = W - Deg;
% 
% K_v = expm(beta * L_v);
% K_v_inv = inv(K_v);
% 
% save('u1_movie_diff_kernel.mat', 'K_v');
% save('u1_movie_diff_kernel_inv.mat', 'K_v_inv');

%*******************************************************************
% Commute Time kernel

% The degree matrix
% Deg = diag(sum(W));
% 
% % The graph laplacian
% L_v = Deg - W;
% 
% K_v = pinv(L_v);
% K_v_inv = inv(K_v);

% save('u1_movie_ct_kernel.mat', 'K_v');
% save('u1_movie_ct_kernel_inv.mat', 'K_v_inv');

%*******************************************************************
% Regularized laplacian kernel

% sigma = 5;
% 
% % The degree matrix
% Deg = diag(sum(W));
% 
% % The graph laplacian
% L_v = Deg - W;
% 
% % Regularized laplacian
% L_v_reg = Deg^(-0.5) * L_v * Deg^(-0.5);
% 
% K_v = inv((eye(N,N) + sigma * L_v_reg));
% K_v_inv = inv(K_v);

% save('u1_movie_regLap_kernel.mat', 'K_v');
% save('u1_movie_regLap_kernel_inv.mat', 'K_v_inv');


%*******************************************************************
% % Diffusion kernel from genre info
% 
% numMovie = size(G, 1);
% numGenre = size(G, 2);
% 
% % The output kernel
% K_v = zeros(numMovie, numMovie);
% 
% % The intermediate symmetric matrix
% W = zeros(numMovie, numMovie);
% 
% % Construct W based on the similarity between each two movies. W(i,j) =
% % sim(i,j) = #common genres between Movie i and j
% 
% for m = 1 : numMovie
%     for i = m : numMovie
%         comGenre = find(G(m,:) == 1 & G(i,:) == 1);
%         W(m,i) = numel(comGenre);    
%     end
%     W(m,:) = W(m,:)/W(m,m); %FIXME: is this normalization necessary?
%     W(:,m) = W(m,:);
% end
% 
% % The degree matrix
% Deg = diag(sum(W));
% 
% % Normalization
% % for m = 1 : numMovie
% %     W(m,:) = W(m,:)/Deg(m,m);
% %     W(:,m) = W(m,:);
% % end
% 
% % Negation of graph laplacian
% L_v = W - Deg;
% 
% % Eigen decomposition of L
% [EigVec, EigVal] = eig(L_v);
% 
% % Exponentialized D
% Deg_exp = zeros(numMovie, numMovie);
% for m = 1 : numMovie
%     Deg_exp(m,m) = exp(beta * EigVal(m,m));
% end
% 
% % Diffusion kernel
% K_v = EigVec * Deg_exp * EigVec';
% K_v_inv = inv(K_v);
% save('movie_diff_kernel.mat', 'K_v');
% save('movie_diff_kernel_inv.mat', 'K_v_inv');